      SUBROUTINE NMSIMP(X,N,FOPT,FUNC)
CNMSIMP     SETS UP WORK AREA AND CALLS NELMED
      IMPLICIT REAL*8 (A-H,O-Z)
      DIMENSION X(N)
      COMMON/BSTAK/NQ,NTOP
      COMMON/BSTACK/A(150000)
      COMMON/BOPT1/NPARS,JH,JFP,JFP1,JSP,JSP2,JA1,JS,JBHHH
      COMMON/BOPT2/ACC,R ,PM1,IVAL,ITERL,ITERC,MX,IER
      COMMON/BOPT3/JSMP,JF,JSPD
      COMMON/BNEL/STP1,VAR,KONVGE,NRST
      EXTERNAL FUNC
      VAR=0.
      IF(KONVGE.LE.0)KONVGE=N+2
C  P(N,N+1), F(N+1)
      NT=NTOP+N*(N+4)+1
      IF(NQ.LT.NT)GOTO 8
      JSMP=NTOP+1
      NTOP=NT
      JA1=JSMP+N*N+N
      JS2=JA1+N
      JF=JS2+N
      NN=N+1
 10   CALL NELMED(N,X,FOPT,A,A(JA1),A(JS2),A(JF),NN,FUNC)
c10   CALL NELMED(N,X,FOPT,A(JSMP),A(JA1),A(JS2),A(JF),NN,FUNC)
      RETURN
 8    IER=-4
      RETURN
      END
      SUBROUTINE NELMED(N,X,FOPT,P,STAR1,STAR2,F,NN,FUNC)
CNELMED              NELDER MEAD(AS47)
      IMPLICIT REAL*8 (A-H,O-Z)
      DIMENSION X(N),P(N,NN),STAR1(N),STAR2(N),F(NN)
      character*4 SR(5),SOURCE
      COMMON/BOPT2/ACC,R ,PM1,IVAL,ITERL,ITERC,MX,IER
      COMMON/BNEL/STP1,VAR,KONVGE,NRST
      COMMON/BPRINT/IPT,NFILE,NDIG,NPUNCH,JPT,MFILE
      COMMON/BSTOP/NVAR1,ISTOP(3)
      COMMON/BPREC/RSMALL,ABSMAL
      EXTERNAL FUNC
C        REFLECTION,EXTENSION AND CONTRACTION COEFFICIENTS.
      DATA SR/4HSTRT,4HRFLC,4HEXTN,4HCTRC,4HALL /
      DATA RCOEFF/1.D0/,ECOEFF/2.D0/,CCOEFF/5.D-1/
      SOURCE=SR(1)
      DEL=STP1
      IRST=0
      JCOUNT=0
      CALL FUNC(X,N,FOPT,*904)
      IVAL=IVAL+1
      IF(IPT.GT.0.OR.JPT.GT.0)CALL OPTOU0(X,N,FOPT,4)
C        CONSTRUCTION OF INITIAL SIMPLEX
 1001 CALL MOVE(X,P(1,NN),N)
      F(NN)=FOPT
      DO 2003 K = 1 , 10
      DO 2 J=1,N
      XJ=X(J)
      X(J)=XJ+DEL
      CALL MOVE(X,P(1,J),N)
      CALL FUNC(X,N,F(J),*2002)
      IVAL=IVAL+1
 2    X(J)=XJ
      GO TO 1002
2002  IVAL = IVAL + 1
      X(J) = XJ
2003  DEL = DEL / 2.
      IF (IPT.GT.0) WRITE (NFILE,2006)
      IF (JPT.GT.0) WRITE (MFILE,2006)
2006  FORMAT ('1INITIAL SIMPLEX COULD NOT BE BUILT WITH 10 REDUCTIONS',
     * ' OF STEP SIZE')
      IER = -12
      RETURN
 1002 IF(ITERC.GE.ITERL)GOTO 900
      IF(IPT.GE.2.OR.JPT.GE.2)CALL OPTOU1(N,ITERC,SOURCE,IVAL,FOPT,X)
      ITERC=ITERC+1
C      LOOP BEGINS HERE
 1000 CONTINUE
C        SIMPLEX CONSTRUCTION COMPLETE
C        FIND HIGHEST AND LOWEST F VALUES. FBAD (=F(IBAD) ) INDICATES
C        THE VERTEX OF THE SIMPLEX TO BE REPLACED.
      IOPT=1
      IBAD=1
      DO 5 I=2,NN
      IF(F(I).LE.F(IOPT))GOTO 4
      IOPT=I
      GOTO 5
 4    IF(F(I).GE.F(IBAD))GOTO 5
      IBAD=I
    5 CONTINUE
      IF(MX.EQ.1)GOTO 51
C   FOR MIN INTERCHANGE IOPT AND IBAD
      I=IBAD
      IBAD=IOPT
      IOPT=I
 51   FOPT=F(IOPT)
      IF(IPT.GT.2.OR.JPT.GT.2)CALL OPTOU1(N,ITERC,SOURCE,IVAL,FOPT,
     * P(1,IOPT))
C        CALCULATE  THE CENTROID OF THE SIMPLEX VERTICES
C        EXCEPTING THAT WITH F VALUE FBAD.
      DO 6 I=1,N
    6 STAR2(I)=0.
      DO 7 J=1,NN
      IF(J.EQ.IBAD)GOTO 7
      DO 75 I=1,N
 75   STAR2(I)=STAR2(I)+P(I,J)
    7 CONTINUE
C        REFLECTION THROUGH THE CENTROID.
      Z=(1.+RCOEFF)/DFLOAT(N)
      DO 8 I=1,N
 8    STAR1(I)=Z*STAR2(I)-RCOEFF*P(I,IBAD)
      CALL FUNC(STAR1,N,CF1,*15)
      IVAL=IVAL+1
      IF(CF1*PM1.LE.FOPT*PM1)GOTO 12
C        SUCCESSFUL REFLECTION, SO EXTENSION
      Z=(1.-ECOEFF)/DFLOAT(N)
      DO 9 I=1,N
 9    STAR2(I)=ECOEFF*STAR1(I)+Z*STAR2(I)
      CALL FUNC(STAR2,N,CF2,*19)
      IVAL=IVAL+1
C        RETAIN EXTENSION OR CONTRACTION.
      IF(CF2*PM1.LE.CF1*PM1)GOTO 19
      SOURCE=SR(3)
 10   CALL MOVE(STAR2,P(1,IBAD),N)
      F(IBAD)=CF2
      GO TO 901
C        NO EXTENSION.
 12   IF(IPT.GT.2.OR.JPT.GT.2)CALL OPTOU1(N,ITERC,SR(2),IVAL,CF1,STAR1)
      IF(CF1*PM1.LE.F(IBAD)*PM1)GOTO 15
      DO 13 I=1,NN
      IF(F(I)*PM1.LT.CF1*PM1 .AND. I.NE.IBAD)GOTO 19
   13 CONTINUE
C        CONTRACTION ON THE REFLECTION SIDE OF THE CENTROID.
      CALL MOVE(STAR1,P(1,IBAD),N)
      F(IBAD)=CF1
C      CONTRACTION ON THE F(IBAD) SIDE OF THE CENTROID.
 15   Z=(1.-CCOEFF)/DFLOAT(N)
      DO 16 I=1,N
 16   STAR2(I)=CCOEFF*P(I,IBAD)+Z*STAR2(I)
      CALL FUNC(STAR2,N,CF2,*17)
      IVAL=IVAL+1
      SOURCE=SR(4)
      IF(CF2*PM1 .GE.F(IBAD)*PM1)GOTO 10
C        CONTRACT WHOLE SIMPLEX
 17   DO 18 I=1,NN
      IF(I.EQ.IOPT)GOTO 18
      DO 181 J=1,N
 181  P(J,I)=.5*(P(J,I)+P(J,IOPT))
      CALL FUNC(P(1,I),N,F(I),*906)
      IVAL=IVAL+1
   18 CONTINUE
      SOURCE=SR(5)
      GO TO 901
C        RETAIN REFLECTION
 19   CALL MOVE(STAR1,P(1,IBAD),N)
      F(IBAD)=CF1
      SOURCE=SR(2)
  901 JCOUNT=JCOUNT+1
      IF(F(IBAD)*PM1.LT.F(IOPT)*PM1)IBAD=IOPT
      CALL MOVE(P(1,IBAD),X,N)
      FOPT=F(IBAD)
      IF(JCOUNT.NE.KONVGE) GOTO 1000
      JCOUNT=0
      Z=0.
      DO 902 I=1,NN
  902 Z=Z+F(I)
      Z=Z/DFLOAT(NN)
      CURMIN=0.D0
      DO 903 I=1,NN
 903  CURMIN=CURMIN+(F(I)-Z)**2
      VAR=CURMIN/DFLOAT(N)
C        VAR IS THE VARIANCE OF THE N+1 FN VALUES AT THE VERTICES.
      IF(VAR.GE.ACC*ACC)GOTO 1002
C        FACTORIAL TESTS TO CHECK THAT FOPT   IS A LOCAL MINIMUM
      DEL=DMAX1(ACC*1.D-3,RSMALL)
      DO 24 I=1,N
      XI=X(I)
      CALL FUNC(X,N,CF1,*23)
      IVAL=IVAL+1
      IF(CF1*PM1.GT.FOPT*PM1)GOTO 25
 23   X(I)=XI-DEL
      CALL FUNC(X,N,CF1,*24)
      IVAL=IVAL+1
      IF(CF1*PM1.GT.FOPT*PM1)GOTO 25
 24   X(I)=XI
C      OPTIMUM REACHED
      IER=4
      IF(IPT.GT.0)WRITE(NFILE,913)
      IF(JPT.GT.0)WRITE(MFILE,913)
      RETURN
C        RESTART PROCEDURE
 25   DEL=DMAX1(ACC,RSMALL)
      FOPT=CF1
      IRST=IRST+1
      IF(IRST.GT.NRST)GOTO 905
      IF(IPT.GT.0)WRITE(NFILE,911)IRST
      IF(JPT.GT.0)WRITE(MFILE,911)IRST
      GOTO 1001
 900  IER=-1
      RETURN
 904  IER=-9
      RETURN
 905  IF(IPT.GT.0)WRITE(NFILE,912)
      IF(JPT.GT.0)WRITE(MFILE,912)
      IER=-10
      RETURN
 906  IF(IPT.GT.0)WRITE(NFILE,910)
      IF(JPT.GT.0)WRITE(MFILE,910)
      IER=-11
      RETURN
 910  FORMAT(27H1ERROR IN SIMPLEX SHRINKAGE)
 911  FORMAT(16H0RESTART NUMBER ,I3,12H STARTS HERE)
 912  FORMAT(40H1NUMBER OF PERMISSIBLE RESTARTS EXCEEDED)
 913  FORMAT(50H1VARIANCE IN SIMPLEX IS LT ACC**2 AND OPT ACHIEVED)
      END
